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FOREWORD 


I 


This report, Volume I of 2 volumes, was jointly prepared by 
the Guidance and Controls Section and the Scientific Programming 
Section of Martin Marietta Aerospace, Denver Division, under Contract 
NAS8-30808. This volume contains the philosophy and the mathematical 
basis of the non-linear programming algorithm underlying the develop- 
ment of the COEBRA program. Volume II is the User's Manual for the 
COEBRA program. The purpose of the contract was to convert the COEBRA 
program from the CDC 6400/6500 digital computer system to the UNIVAC 
1108 at the George C. Marshall Space Flight Center and to provide a manual 


and instruction on the use of the program. 


This contract was performed 


from September 1974 to August 1975, and was administered by the National 
Aeronautics and Space Administration. George C. Marshall Space Flight 
Center, Huntsville, Alabama, under the direction of Mr. D. K. Mowery, 
Dynamics and Control Division, Aeroastrodynamics Laboratory. 
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CHAPTER 1 


THE DESIGN PROBLEM 


1,1 Introduction 

The detailed design of the autopilot for a large elastic booster is 
at best a tedious procedure. The function of the autopilot is to insure 
satisfactory flight from liftoff to orbital entry or touchdown, depending on 
the specific mission involved. In general, the autopilot works in conjunction 
with, but independent of, the booster guidance system. It is primarily a 
stability controller functioning to stabilize the vehicle in order that (1) 
the vehicle response to guidance commands is not compromised, (2) trajectory 
dispersions resulting from both internal and external disturbances be minimal, 

(3) structural integrity of the vehicle be maintained. 

Typically, the autopilot designer is provided with a mission profile 
including nominal trajectory data and vehicle/payload dynamics. A critical 
element of the trajectory data will include specific flight conditions that must 
be met at specified points along the trajectory. These points are usually points 
at which a change in autopilot parameters, filter gains and/or compensation time 
constants is required, or at which either critical vehicle dynamics are encounters 
or trajectory maneuvers are scheduled to occur. They are identified in this 

volume as trajectory events or trajectory time points. 

The control conditions to be met at these critical trajectory points 
can be specified in terms of stability criteria based on an open loop frequency 
response analysis or a closed loop time response and are often specified as a 

mixture of both. 

The common stability "performance measures" are defined in either the 
real frequency domain or in the time domain and are derived from the basic 
properties of the governing differential equations or the transfer functions 


(Laplace transforms of the input-output differential equations) for the 
system under consideration. 

1^2 Analysis in the Frequency Domain 

Without going into great detail, it is assumed that the reader has 
covered the ground before, recall that the general solution of a linear 
differential equation(s) can be written as an exponential function or a 
linear combination of exponential functions with exponents >it. The >i can 
be real, complex, or imaginary, depending on the nature of the specific equation(s) 
being solved. If the real part of any one or more of the ^i is positive, the 
system governed by the differential equation is said to be unstable. If the 
real part of any one or more of the ^i is identically zero, the system is said 
to be oscillatory. If the real part of all the >i are negative, the system is 
said to be stable. 

The numerical values of the ^i are derived from the characteristic 
equation of the system, and are the poles of the system transfer function. We 
are thus led to a consideration of the location of the system poles and a 
characterization thereof that will yield the greatest amount of information 
regarding the overall system behavior. The hyquist stability criteria provides 
such a characterization and is typically used in the analysis of space vehicle 

autopilots . 

Let us consider a simple control system as a preliminary. Figure 1-1 
exhibits such a system. In this system we relate to the control elements, 
to the controlled system or "plant," H to the feedback elements, and and 
to the forward and feedback path gains. To simplify what follows, we let 
G = Gi G 2 and the closed loop transfer function becomes 

^ G(s) 1-1 

R " 1 + Ki K^G(s) H (s) 


2 


Figure 1-1 Typical Control System 


In the application of Nyquist criteria, we treat the complex variable 
3 = o- (li as a pure imaginary, or <>'=0 and examine the characteristic equation 

1 -I- G(s) H(s) = 0 ^'2 

Rather than treat the real frequency behavior of the closed loop system, we 
treat the open loop function, G H, relative to the complex quantity, 

-1=1 e^^. The Nyquist criteria states that, if in the complex frequency 

domain, we map 

F(s) = G(s) H(s) . , . . - c , . o . . . . ^-3 

as s traverses the contour of Figure l-2a in the s-plane the mapping of F(s) does 
not encircle or pass through the -1 point in the F(s) plane, none of the roots 
of F(s) lie in the right half s-plane or along the imaginary axis of this plane. 
Thus, in the time domain, the real part of all the^i are all less than zero and 
the system is stable. 

Typically, along the contour T ^ of Figure 2a, the quantity F(s) vanishes 
and we are left only with F(s) mapped as s traverses and F^, usually in 

that order. We further observe the generally F(j<fl) for positiveuis the mirror 










image of F(-ja>) and thus we need only map F(jtt)) for the real frequency 

domain. Thus, wc can actually observe not only the absolute stability but 
also the relative stability through examination of the teal frequency response. 
Consider F(j<^) = Re (K^K^GH) -f- j (K^K^GH) 1-4 

The relative stability can be exhibited by the gain margin, 




KiK^GH 


-1 


- 1 


1-5* 


where is that frequency for which the argument of F(ju)) or the phase angle 
of Kj^K^GH is 7T radians or 180 degrees, and the phase margin 


p- arg F(j(i))| radians 
(*> « 


l-6a 


= 180 - arg F(jw) degrees , . l-6b 

I * "c 

where to > the crossover frequency, is that frequency at wiiich F(j<o)l = 1, 
c I I 

Additional useful information regarding equivalent damping and the location of 

r45i** 

the least stable poles of F(jto) are discussed in the open literature.*- 

Figure 1-3 illustrates these margins in the conventional manner. 

The advantage of the use of the Nyquist diagram in the real frequency 
domain in the design of elastic booster autopilots lies in the fact that the 
rigid body, aerodynamic and flexible body bending mode frequency responses are 
relatively well separated and thus lend themselves to definitions of separate 
gain and phase margins, as shown in Figure 1-4, an unfolded Nyquist diagram with 
frequency as a running parameter and magnitude and phase as coordinates. 


*Gain margin is defined as the gain increase required to drive the system into 
*instability or K |F(j(o^)| = 1 and K' = 1 f(J(^) | , 

Thus, in a vector sense, the margin is given by 1-5, 

**figures in brackets refer to references in the Bibliography. 
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Imaginary 



\ arg FCjtu) 


= 180-arg F(ju>)| 


Figure 1-3 Definition of and in terms of Nyquist plot 


1.3 Analysis in the Time Domain 

In the time domain the performance measures are usually defined in 

terms of the system response to a step function Input. The generally accepted 

performance or stability measures, based on a closed loop system, are the rise 

time, T^, settling time, T^, time to peak value, time delay, T^, percent 

overshoot and steady-state error. The rise time is commonly defined as the 

time required for the system to go from 0.1 to 0.9 of the final value. The 

time to peak, T^, is the time to reach the first peak of overshoot and represents 

one half cycle of the frequency of the damped oscillation typical of systems. 

The settling time, T , represents the longest exponential decay of the system 

s 

and is defined as the time required for the system to settle within a band 0.95 
to 1.05 of the desired final value. The delay time, T^, is the time required for 
the system to reach 0,5 of the final value. The overshoot is the amount the first 
half cycle of the damped oscillation exceeds the desired final value and the 
steady state error is the difference between the final value of the system 
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FIGURE 1-4 Open Loop Gain/Phase Plot - Unfolded Nyquist Diagram in 
Rectangular or Nichols Chart Coordinates 




response and the desired value (typically 1). Figure 1-5 illustrate.y these 
parameters . 


It is of interest to note that many high order systems can be approxi 

mated by an equivalent second order system* In any case, the time response 

is governed by what are kno\^m as dominant poles. To better appreciate what 

is meant by a dominant pole, consider a closed loop system function of the 

general form 
J 

K TT + 

c(s)= 1-7 

le,.i J . « 


where 2i -h + 2 p + n, or simply there is an excess of poles over zeros. 

This enables us to form a partial fraction expansion of C(s) in elementary 
terms of the form 

L^l K-.\ ^ ^ ^ 




This, in turns, inverts to: 


J = ' 






x:= ‘ 




Observe now the coefficients of t in both the exponential and cosine terms. 
These are derived from or are identical with the system poles with a step 
function, M (t) , input. Those terms of equation 1-9 that give rise to the 


significant portions of c(t) stem from the dominant poles. That is, those 
terms dominate the performance of the system in the time domain. Their locations 


are clearly functions of the system parameters and are closely related to the 

[^ 6 ] 


gain margin, phase margin, rise time, settling time, etc,^ The importance of 
these poles is discussed further in Appendix B in conjunction with examples. 


1 ,4 The Autopilot Problem 

The foregoing indicates that the design of a control system is a relatively 
s Lraight- forward engineering problem. However, when the numerous inputs and/or 
disturbances, multiple loops and/or multiple control effectors, or a combination 
of these are included, as in a typical flight vehicle autopilot, the problem becomes 
both increasingly difficult and tedious, with results somewhat less than 
immediately visible, A typical autopilot may have as many as 6 distinct feedback 
paths with up to 50 constants that must be established in a reliable and economical 
manner. A mechanization of the iterative design procedure becomes highly 
advantageous . 

Consider a typical autopilot that is to be both gain and phase stabilized. 
The vehicle flight plan is designed to minimize aerodynamic loading. The vehicle 
is launched vertically to eliminate or minimize (lateral stabilizing load require- 
ments)on the autopilot. Shortly after lift-off the vehicle is rotated (pitched 
over) to assume a zero- lift (zero angle of attack) flight condition. The pitch 
rate is then maintained at the flight path turn rate to obtain a. zero angle of 
attack. 

The autopilot designer is supplied with this nominal trajectory as well 
as nominal vehicle dynamics including rigid body, bending mode, and fuel slosh 


10 






















characteristics. Assuming that a standard autopilot configuration has 
already been established, the usual case today, the first step is to 
initialize the parameters used within this standard configuration. The initial 
gains for compensation are derived from the thrust to inertia plot for the 
specific vehicle, based on a gain margin of about 6 db. The initial time 
constants are derived from the structural bending mode frequencies. 

Figure 1-6 illustrates the commonly accepted configuration used. The nominal 
trajectory and vehicle dynamics (mass and inertia properties) are generally 
available as the results of a computer run and thus can be provided as a deck 
of punched cards, magnetic tape or other computer compatible format. Using 
these data, and the standard configuration with the initial gain and time 
constants, the frequency response is plotted. A standard computer routine 
based on a block diagram input format according to the block diagram of Figure 1-7 
is used where each block represents an elementary transfer function of the 


general form 


T-(5)= 

^ 






±-io 




It is to be noted here that T(s) of equation 1-10 is identically the form 
of G(s) H(s) of the total open loop function of equation 1-3 where for the 
elementary blocks, m, 1, j, n, q and v are of lower magnitude than for the 
total open loop function. The resulting frequency response is typical of that 
of Figure 1-4. The response curve is examined for gain and phase margins. The 
designer then selects new gain constants and/or time constants and repeats the 
computer run until a satisfactory response satisfying all gain and phase margin 
constraints is obtained. 
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1,5 A Simple Example 

As an example. «e treat the system o£ Figure 1-8 with a specification 
that the velocity constant, K^, be 100, the crossover frequency. w^>10 rad/sec, 
the gain margin. >10 db and the phase margin, M^>45°. A further constraint 

is that the system be a unity feedback configuration. 



160 


> 

s (S-F8) Cs +20) ! 



Figure 1-8 Plant of Example 

Recall that the velocity constant, K^, is defined as Che ratio of input 
rate, r' , to the steady state error, E, and it can be shown that 

K = lim s = lim s G(s) 

V s _^0 ^ s -^0 

Plotting the conventional Bode diagram and the commonly used gain-phase diagram 
as in Figures 1-9 and 1-10, we observe that the 3rd order plant, in a unity 
feedback configuration is unstable. The gain crossover frequency, is about 

21 radians per second, at which frequency the gain, K, is 1 and the phase is 
approximately 210 degrees. The Nyquist diagram, of which Figure 1-10 represents 

a small segment, has the general shape shown in Figure 1-U. 

The system is clearly unstable and series compensation in the forward 
loop is indicated. After some consideration, a lag-lead compensation network 

the following configuration (Figure 1-12). 


is arrived at, resulting in 
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Figures 1-13 and 1-14 show the Bode and Gain-Phase (Nyquist) plots 
for the compensated system. The phase margin is about 54° and the gain 
margin is about 15 db. Clearly, the advantage of a computerized solution 
to this simple system is doubtful. However, the number of parameters, 
feedback loops and variables in the autopilot problem evidences this advantage 
to a marked degree. Recall also that only the starting point and a usable end 
point were demonstrated for this simple example. Furthermore, we observe that 
the solution is not optimal in that the gain margin exceeds requirements by 
about 5 db and the phase margin is about 9° more than required. 

We are thus led to develop a computer program that will not only 
produce the system frequency response in an efficient manner but also will 
enable an automatic search for near optimum gains and time constants based 
on a gain margin and/or phase margin constraint or alternatively, the time 
response and its related performance measures. 
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Figure 1-14 Gain-Phase Plot as a Function of Frequency for Bode Plot of Fig. 1-13 
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CHAPTER 2 


THE DESIGN ALGORITHMS 


2.1 Linear Progranrnilng 

With gain and phase margins, or any other performance measure, as 
the autopilot design criteria, the selection of filter time constants and 
gain factors can be realized through the use of the linear programming 
technique. With this thought in mind, we examine briefly the linear pro- 
gramming approach before we apply it directly to the autopilot problem. 

We treat the linear programming technique in two dimensions in order 
to clarify the approach. Linear programming addresses the solution of a set 
of m linear equations in n unknowns subject to constraints on the unknowns and 
a minimizing functional. Consider 

A X = 2-1 

where A is an m x n coefficient matrix, X an n vector of unknowns and B is an 
n vector of constants. We impose the following constraints: 


X j ^ 0, j==l,2, 


n 


2-2 


and 


n 

£ 


j=i 


c . X , = mintmum 
J 3 


2-3 


The details of the linear programming method are set forth in the open literature 
and will not be treated fully here; rather, an example will serve our purposes. 
Given the set of inequalities: 


x^ + 3x^ 

<10 



X. + x„ 

<6 1 



1. 2 

\ • • « • • < 

• • • • • fl 

2-4 

Xi - X2 

^2 ( 



x^ -h 

>6 ' 




23 






First we alter 


1 

with and non-negative, minimize the function -x^ - 2 X 2 - 
our minimal function constraint to maximize the function x^+ 2x2* The first 
step is to bound the required solution. This is done by changing all the In- 
equalities to equalities and plotting each equation of the set 2-4, recognizing 
that the equalities represent either upper or lower bounds on the inequalities. 

The region below the lines -x^ + 3 X 2 = 10 and x^ + X 2 = 6 contain admissable 
values of x^ &. X 2 « The entire first quadrant satisfies the requirement that 
x^ and X 2 be non-negative. The region above the lines x^ - X 2 =2 and x^ + 3 X 2 = 6 
may contain admissable values. Thus, the region bounded by these lines contains 
an admissable value of both x^^ and X 2 as indicated in Figure 2-1. Overplotting 
the function to be maximized we write 

+ 2 X 2 = Constant 2-5 

and plot a series of lines overlaying the region of admissable solution to 
obtain the admissable solution to the set 2-4 that maximizes the function ‘*'^^2* 
Inspection of figure 2-2 indicates that the macimum admissable value of x^ + 2 X 2 
is 10, from which we see the desired solution to the set 2-4 is = 2, X 2 ” 4. 

This simple two-dimensional (two variable) problem is almost trivial. 
However, as the dimension n increase::, the simple 2-space geometry can no longer 
be used. The concept, for an n-dimensional space, is still the same. Using a 
notational scheme that is not dimensionally limited and compacts the statement 
of the linear programming problem, we write: Determine the column vector X, a 


solution of 

A X = B 2-5 

which maximizes 

Y = X 2-6 

where X^O 2-7 


We denote 2-6, where C is the transpose of the column vector C, the objective 
function and 2-5 the constraints. The non-negative property, 2-7, is fundamental 
to the linear programming algorithm, 

24 















How do we apply this procedure to the autopilot problem? We consider 


thvi gain margin, phase margin, rise time, settling time, maximum overshoot, etCp 
as performance measures or objective functions, and constrain their values to 
an acceptable range. These acceptable values usually will err on the conserva- 
tive side in order that unforeseen parameter variations will not allow excursions 
outside of "safe" values for the performance measures. 

In the frequency domain, it is a relatively simple matter to identify 
the frequencies at which the performance measures are obtained. For example, a 
gain margin exists at those frequencies for which the phase of the system is 180 
degrees, and a phase margin exists at frequencies for which the gain is zero 
decibels (See Figure 1-4). The frequencies at which pertinent performance measures 
are found are generally well separated, and in sequential order with increasing 
frequency, e.;abling a simple search routine to identify each measure in sequence. 

In the time domain, a similar search technique permits identification of 
rise time, settling time, overshoot and Che maximum deviation from final value,* 

We thus treat these performance measures here in a general way to 
illustrate the design technique that has evolved into the COEBRA design algorithms. 

Let the vector X represent that sot of autopilot parameters whose values: 

1) Determine directly or indirectly the value of one or more |j 

■? ' 

performance measures, 

2) Can be varied in order to achieve an "optimum" performance measure 
in some reasonable sense, 

3) Can be changed, periodically, during flight to maintain satisfactory 
vehicle stability, 

* Due to structural reasonances, it is possible to have a time response in which 
the maximum excursion and the first half cycle overshoot are not coincident 
in time. 


X 




as a function of this parameter 


Thus we identify each performance measure 
set, assign a nominal value to them, and use these values to initialize an 
autopilot designo A frequency response or time response is then calculated 
and the set of performance measures examined. Assuming that a specific 
performance measure specification is not made, a single element of the parameter 
vector is changed and the frequency response or time response is repeated. The 
method of finite difference is then utilized to approximate the partial derivative 
of the performance measures with respect to this parameter to implement the 


performance measure constraint equation 

n 


M 


■ f ^ -f- A x> ‘ ax, t 


6-^2; • " V 


yn 




where M. is the i-th performance measure, the initial value of the parameter 

vector, M, the specified performance measure minimum and Ax^ ^ ‘ ^jo* 

is 

Rewriting 2-8 with the substitution for Axj 






Irisov;'-' 



2-9 


. . yn 


it is possible to invoke more than one mission flight condition (trajectory 
time event) in evaluating the specific performance measures, M., under con- 
sideration. 

A fundamental restriction on this design algorithm is that the user 
must specify: 1) the autopilot configuration within the limits identified by 


*2-8 is actually a truncated Taylor series expansion of the performance 
measure, M. , about its initial value for each interation cycle. 



the block diagram of Figure 1-7 or in a particular matrix format (see 
Appendix D) ; 2) the number of gains and filter time constants (elements of 
the vector X) which are to be treated as variables and; 3) those autopilot 
parameters which are to be treated as constants for each iteration cycle. 

The algorithm then operates on an individual element of the vector X under the 
constraint 


X..-<i X. < X. , 

jl jn 


j = 1, 2. 


2-10 


where the individual may include several different flight conditions 
(trajectory times or vehicle states) or any one or more x^ may remain constant 
over several flight conditions or trajectory time events. 

2,3 The Objective Function 

We now formulate the objective or cost function in such a way that as 
the performance measures improve for one flight condition, the trend is for 
them to improve for all flight conditions, and at the same time each structural 
bending mode resonance will be forced to occur in the neighborhood of zero 
degrees phase. This is done by forcing the partial derivatives in the constraint 
equation to increase. In the format of equation 2-6, we thus write 


'l = •? W^ (i, t) • (i, t) • S (i) • (i» t) . Xj 2- 


10 


where : 


1) is the summation over all autopilot variables 

2) X is the summation over all flight conditions (trajectory time 
points or vehicle states), 

3) ^ is the summation over all performance measures at all flight 
conditions. 

4) (i, t) refers to a weighting factor. For each performance measure, it is 
simply a ratio of the desired measure over the actual measure. 

Hence, if a performance measure requirement is not met, (i, t) 
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4) (Continued) 


will be greater than unity. It becomes less than unity when a 
margin exceeds its desired objective. It is noted at this time 
that in the expression for yj^, i also indexes the phase angle at 
which each structural bending mode reasonates. For these values 
of i, the partial derivative indicates the rate of change of 
each modal peak phase with respect to each autopilot variable, 
and (i, t) is written so that the algorithm will attempt to 
force each mode to resonate near zero degrees phase, (i, t), 

will be large for modes that resonate near 180 degrees, and zero 
for modes that resonate at zero degrees. For some arbitrary angle 
like 90 degrees, (i, t) can equal unity. 

5) (i, t) refers to a weighting factor that might be selected by 
the user. This would give the user the capability to eliminate 
certain performance measures from the optimization process or to 
emphasize other measures.* 

6) S(i) refers to a scale factor. It serves to scale the margins and 
modal peak phases so that phase margins and gain margins can be 
optimized together. For example, it might be desired to equate a 
five degree increase in the rigid-body phase margin with a one 
decibel (12,27.) increase in the rigid body gain margin. Similarly, 
in the time domain, a correlation between rise time and overshoot 
may require scaling. 

In summary, is a "weighted" linear combination of the "positive" 
changes in each performance measure. Note that this linear combination can 


*W^ (i,t) and (i,t) arc initially preprogrammed, After an initial run of 

COEBRA, the user can, at his discretion, input new values for these weightin; 
functions. 


30 










I 
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incorporate performance measures from all of the vehicle states (flight 
conditions) that are being designed together. The design algorithm will 
maximize ^and hence seek to maximize all stability margins and/or time 
domain performance measures, and seek to force all modes to resonate near 
zero degrees phase) in the oresence of the constraint matrix which includes 
constraints on each individual margin and each autopilot variable at each 
time point. The advantage of this algorithm lies in the fact that the constraint 
equations can specify the minimum requirements for each performance measure while 
the cost function seeks to maximize the performance measure, 

2 Load Relief Cost Function 

Structural bending moment loads on a launch vehicle are largely due to 
axial acceleration, aerodynamic loading, and control device deflections [Harris, 
15], Obviously, the booster autopilot can do little to affect axial acceleration, 
and therefore the main objective of a so-called load relief autopilot is to 
reduce aerodynamic loading due to angle of attack and to keep control device 
deflections to a minimum. 

Hence, for this design algorithm, when the objective is to maximize 
structural bending moment load relief capability, the cost function is comprised 
of the response of the angle of attack (^) and the control deflections ( 5 ) due to 
the wind dorcing function (^w) • When the cost function is maximized, the peak 
values of ^ and 5 are minimized, 

A separate transient response routine is used to calculate the peak 

value of angle of attack (fl ) and control deflection (6 ) due to . As with 

P p 

stability margins, the method of finite differences is used to compute the 
first partial derivatives of and with respect to the autopilot variables. 

The cost function is then formed from the first order terms of the Taylor Series 
expansions of and 5 ^ about their nominal values. 
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As with the stability margin cost function, the lodd relief cost 
function (y 2 ) is a weighted linear combination of the variable portion of 
the first order terms in the Taylor Series, or 


z z 


(t) * 


1£e 


(t) 


ax. 

j 


+ «2 * 


ax, 


* X 


2-12 


In the above expression: 

(1) 'prefers to the summation over all the autopilot variables; 

(2) refers to the summation over all the vehicle states; 

(3) (t) and refer to weighting factors that may be 

specified by the user. 

When maximizing load relief capability, the design algorithm will 
maximize the negative of ii^ the presence of the constraint equations on the 
minimum allowed gain/phase stability margins and on the allowed ranges of the 
individual autopilot variables. Note that multiple time point design is handled 
just as it is when maximizing other performance measures. Some final notes on 
the load relief cost function are now listed. 

Since the so-called "rigid-body" (as opposed to flexible-body) angle 
of attack (/3) and control deflection (8) arc the principal factors in determining 
structural bending moment loads, it is felt that only the rigid-body airframe 
equations of motion [Harris, IS] need to be used in the transient response routine 
that is used to calculate angle of attack and control deflection. Note also that 
these rigid-body airframe equations can include planar coupling (e.g., between 
the yaw and the roll planes), and hence the cost function can include control 
deflections from several planes (e.g., the yaw plane control deflections (8^) and 
the roil plane control deflections. (8^) ). 
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The wind forcing function can be a series of steps and/or ramps 
that approximate the commonly used synthetic wind profile [l5]. The wind 
forcing function can also be stochastic, and the design algorithm then 
minimizes the rms values of )3and s- This is done via Wiener's theorem and 
the filtering property of power spectral density functions [see Appendix c]* 

2„5 General Flow Chart 

Figure 2-3 is a general flow chart summarizing the main steps involved 
in the algorithm. It shows the general flow from the initial autopilot for 
each iteration through the following routines: 

(1) The routine that generates the frequency response and finds 

'the stability margins, and the routine that generates the transient 
response and finds peak ^ and S, rise time, settling time, overshoot 

etCo 

(2) The routine that computes partial derivatives; 

(3) The routines that set up the linear programming problem and solve 
it; and 

(4) The routines that determine whether the design is complete. 

If the design is not complete, another major iteration is begun with 
the best answer obtained in the previous iteration. In other words, the problem 
is relinearized about the best answer of the previous major loop, and another 
cycle through the major loop is performed. This iterative process continues 
until the local optimum is found. 

Note that this design process satisfies the five main elements of the 
design criteria 

(1) The method directly treats stability margin requirements and 

objectives, and structural bending moment load reduction and/or 
the time domain performance measures. 
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(2) 


I 



i 


The method directly handles the user-selected autopilot 
configuration. 

(3) The method directly handles the multiple time point (flight 
condition) design pr ’em. 

(4) The method is not limited by the order of the system. Note 
that since this is a parameter optimization routine, the order 
of the autopilot does not necessarily increase with an increase 
in the order of the fixed parts of the system. 

(5) The method can design either a digital autopilot or an analog 
autopilot.* 


2.6 S election of Iterative Step Size 

Referring now to equation 2-3 and the sample problem, equation 2-5, 
we observe that to establish an optimum admissable solution, the function to 
be maximized was allowed to take on successive ’’constant" values. The speed with 
which the algorithm converges to a solution is dependent on the magnitude of 
the change in this constant value as the iterative process is carried out. 
Re-casting equation 2-10, to reflect the step-change in each autopilot variable, 

we obtain 


MAX ^1 + p) 
j = 1, 2, . 


-1 




jo jmin 




J 


n 


2-13 


where : 


X and X refer to the minimum and maximum values ever 

^ ' jmin jmax 

allowed for X,. 

J 


*Under this submittal, the formal digital algorithm is not included.^ However, if 
the bi-linear transformation is applied to the pulse transfer function, the 
resulting w-plane description is directly analogous to the s-plane description 
contained herein. Thus, with the w-plane equations, COEBRA can be used directly 
to solve the digital autopilot problem. 


35 




( 2 ) 


1 

t 


X refers to the initial value of X on each iteration. Note 
jo J 

that X is the point about which the partial derivatives are 

jo 

computed, and about which the Taylor Series is expanded, 

(3) P refers to the autopilot variable step-size for each iteration. 

In words, if and X.^^ are not encountered on a particular 

iteration, the above constraint equation says that X^ is allowed to vary no 
more than about + P% from X.^ on any iteration. Since it is desirable to 
maximise the step size on each iteration, thereby getting the maximum "mileage" 
out of each set of partial derivatives, it is desirable to have a Minor Loop that 
increases the size of P until improvement in that "search direction" is no longer 
possible. In other words, Che Minor Loop serves to maximize the autopilot 
variable step-size. In maximizing P, the Minor Loop uses two "indicators": 

(1) a counter that keeps track of the number of performance measures that are 
already met, and (2) a f igure-of-merit that is a linear combination of their 
actual values. If the number of "met measures" increases, obviously the value 
of P can be increased. If the number of "met measures" does not change, the 
flgure-of-merit is used to decide whether P can be further increased. In other 
words, the measure counter is used to reward those steps that result in an 
increase in the number of "met measures." Conversely, the counter prohibits 
those steps that result in a loss in the number of "met measures." Finally, 
the figure-of-merit is used to break ties when the measure counter does not 

change from one step to another. 

The Minor Loop serves to either keep the problem linear on each major 
iteration, or to take advantage of the neglected nonlinearities when they 
might be helpful. In other words, the Minor Loop serves to keep the nonliuearities 

from "hurting" the steady convergence to a local optimum. 

A major benefit of the Minor loop is that it allows the algorithm to 
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converge steadily to an "interior" optimum. 

Since the solution to the linear programming problem always lies at 
a vertex of the admissable region defined by the constraint equations, it is 
the Minor Loop that allows the algorithm to converge to a local optimum that 
is interior to the performance measures constraint equations. 

^ ^ The Inner Loop 

The second part of the step-size optimization routine can be illustrated 
by the following detailed expansion of a particular performance measure constraint 
equation 

* fX. - X.^) (i = 1, . . m) 

2-14 

In the above expression, there are two cases for M. , 

IS 

(1) If the particular pertormance measure is already met, then for 
the next Iteration, 

== SPEC(i) 

where SPEC(i) is the minimum allowed value for the i*" performance 
measure 

(2) If the particular performance measure is not yet met, 

M. = M (X ) + STEP * |SPEC (i) - M (X )] 

IS io o L io o J 

Before defining the purpose of the equations for M. , note that in the 

X s 

second equation, (a) if STEP = 1, M, = SPEC(i), and (b) if STEP = 0, M, s= 

X S X s 

M. (X ). 
lO ^ o' 

For a given value of P (autopilot variable step-size), there may not 
be a feasible solution to the linear programming problem if the present autopilot 
does not meet all of the performance measure constraints. In other words, 
the feasible region defined by the measure constraint equations may not overlap 


M. (X ) 
10 o 


+ L 



the feasible region defined by P. By automatically reducing the value of 
STEP, the measure constraints are "loosened,” until a feasible solution is 
possible for a given value of P. In this way, the value of P can be increased 
in a steady and rational manner, and the algorithm will be allowed to converge 
to a solution in a progressively improving manner, 

2.8 Graphical Illustration of Step-Size Optimization 

Figures 2-4 through 2-9 graphically illustrate the mechanics and 
the interaction of the two step-size optimization routines (the Minor Loop and 
the Inner Loop), 

Figure 2-4 shows a hypothetical two dimensional condition that might 
exist for a rigid-body autopilot design problem. Figure 2-4 is a plot of the 
attitude error gain (K^) versus the attitude rate gain Plotted on the 

figure are three nonlinear stability margin constraint equations: (1) the aero- 

dynamic gain margin; (2) the rigid-body phase margin; and, (3) the rigid-body 
gain margin. Figure 2-4 also shows where the "true” local optimum condition 
might be, where the objective is to maximize stability margins, and where all 
three stability margins are equally weighted. Obviously, the "true" optimum 
for this hypothetical case lies inside the feasible region where all three 
margin requirements are satisfied. The figure also shows what might be the 
"first guess" or initial condition on and 

Figure 2-5 shows what the constraints might look like when they are 
linea-izet' about the initial condition. The figure also shows the linearized 
interior optimum, where again, all margins have been equally weighted. Note 
that the linearized optimum is not the same as the nonlinear or "true" optimum 
for this inir.ial condition. Finally, Figure 2-5 shows the slope of the linearized 
cost function (y) and the direction in which it increases. 
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Figure 2-6 illustrates the feasible region defined by the autopilot 
variable constraint equations on and for step-size #1 (denoted 
This feasible region does not overlap the feasible region defined by the 
stability margin constraint equations. With the initial condition, the rigid 
body phase and gain margin constraints are satisfied, but the aerodynamic gain 
margin is not. Hence, the design method enters the inner loop, and relxes the 
aerodynamic gain margin constraint until a feasible region exists for both the 
margin constraints and the autopilot variable constraints. This relaxation is 
accomplished by reducing the parameter denoted as STEP. 

When STEP is unity, no relaxation exists. When STEP is reduced to 0,8, 
the aerodynamic gain margin constraint is relaxed enough so that a feasible 
region exists. When STEP is 0,8, this means that an "80% improvement" is 
required for the margin that is not yet satisfied. This so-called required 
margin improvement" becomes very important when the optimum is exterior to ti.e 
feasible region. This is the case most of the time for launch vehicle autopilot 
design. 

As indicated on Figure 2-6, the optimum solution for the first step of 
the minor loop exists at y^. Comparing y^^ with the nonlihear cost function and 
constraints shown on Figure 2-4, it is seen that stability at y^ is better than 
at the initial condition. This ’’improved stability" is indicated by the figure- 
of-merit which, as discussed previously, is a linear combination of the stability 
measures. Note that the so-called "met-measure counter" indicates that at y^, thei 
are still only two margins that are satisfied. Note that the measure counter and 
the f igure-of-merit are formed from an actual evaluation of the frequency response 
In other words, they are not computed from the linearized cost function and the 
linearized constraint equations. 

Hence, since y. is better than the initial condition, the design process 



advances, using the same set of partial derivatives, and hence the same 
linearized cost function and margin constraint equations that were calculated 

at the initial condition. 

As shown in Figure 2-7, the design method now increases P from to P^. 
Figure 2-7 shows the feasible region defined by the autopilot variable constraint 
equations for P^. An overlap exists between the feasible regions defined by the 
margin and autopilot constraint equations, and hence the inner loop need not be 
used. For P 2 the optimum solution exists at y^. By comparing ^ with the non- 
linear constraint equations of Figure 2-4, the margin counter indicates that 
there are now three margins that are satisfied. Since improved stability has 
again been achieved, P is further increased from the original initial condition. 

As shov/n in Figure 2-8, P is now increased to P^. For this step, overlap 
also exists, and the optimum solution is at 73 . By comparing 73 to the nonlinear 
constraints of Figure 2-4, it is seen that the rigid-body phase margin require- 
ment is no longer satisfied. The margin counter indicates that only two margins 
are satisfied by 73 . Hence, 73 is not as good as y^ and P must be reduced. 

Figure 2-9 shows the autopilot constraint equations for P^, where P^ < 
P^.^. P 3 . Again, the inner loop is not needed, and the optimum solution exists 
at y^. By comparing y^ to Figure 2-4, the margin counter shows that there are 
three margins that are met. But, the f igure-of-merit shows that 73 is better 
than y^. Postulating that the difference between P^ and P^ is less than some 
convergence criterion, the algorithm stops this so-called major iteration at y^. 
The values of and at y^ become the initial condition for the next major 
iteration. At the problem is relinearized. A new set of partial derivatives 

is computed, and a new cost function and new constraint equations are formed. 

As the design progresses, the linearized optimum gets closer and closer to the 
nonlinear or "true" optimum. It will be observed in the following that the 


43 


convergence criteria can be used to terminate this iterative design process. 

Table 2-1 summarizes the results of Figures 2-4 through 2-9, These 
figures have been used to demonstrate steady convergence to a local interior 
optimum. We will next consider convergence to a local exterior optimum. 

As a final note, Figures 2-4 through 2-9 demonstrate that this algorithm 
does not require that the initial condition lie within the feasible region, 

2,9 Convergence to an Exterior Optimum 

We now illustrate how the algorithm converges to an exterior or "con- 
strained*' optimum. Figure 2-10 shows a case that might exist when optimizing 
load relief capability since for this phase of design, the optimum solution 
almost always is exterior to the feasible region defined by the margin constraint 
equations. Figure 2-10 is a hypothetical two-dimensional case where (1) the non- 
linear margin constraint might represent the so-called aerodynamic gain margin, 

(2) might represent the attitude error gain, and (3) X^ might represent the 
so-called load relief loop gain. Figure 2-10 also shows the nonlinear constrained 

optimum. 

In Figure 2-10 the initial condition on X^ and X^ is outside the feasible 
region. Figure 2-11 shows what the margin constraint and the nonlinear optimum 
might look like when the problem is linearized about the initial conditioi. 

The first step of the algorithm is to "get feasible" and Figure 2-11 ^ j- 11 show 
that in so doing, the algorithm still attempts to approach the optimum. 

Referring to Figure 2-11, after a series of iterations through the 
minor and the inner loops, the solution is shown to exist at 7^. With the 
linearized margin and cost function as shown, this is the best solution this 
major iteration can achieve without violating the nonlinear margin constraint. 

At y^, the problem is relinearized as shown in Figure 2-12, This 
figure shows that any step in the direction of the gradient to the linearized 





Table 2-1 Summary of Figures 2-4 Through 2-9 
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Y serves as initial condition for next major iteration since F^-P^ satisfies convergence criteria 










optimum would yield an unfeasible solution. Referring to Figure 2-12, again 
after a series of iterations through the minor and the inner loops, the 
solution exists at With the linearized margin and cost function as shown, 

this is the best solution this major iteration can achieve without violating 
the nonlinear margin constraint. With each major iteration, the linearized 
constrained optimum approaches the nonlinear constrained opuimum. 

We have illustrated how the algorithm first "gets feasible," and then 
moves along or parallel to a constraint for the case of a constrained optimum. 

Example 6 of Appendix B will dramatically demonstrate this situation. 

Because of the weighting factors in the cost function and figure-of-merit, 
situations with an exterior optimum can also exist when optimizing performance 

measures . 

2 . 10 Termination 

We now discuss the two ways in which the design process can be terminated. 
The first way might be referred to as self-termination, where the algorithm 
finds a local optimum and can achieve no improvement over the initial autopilot 
for a given iteration. The so-called measure counter and the £igure-of-merit 
define this "improvement." Recall that the counter "counts" the number of "met 
measures," and the figure-of-merit is a linear combination of the actual performance 
measures. The counter never allows the algorithm to lose a measure, and an 
iteration is considered better if the counter increases. The figure-of-merit 
is used to break ties in the counter. An iteration is considered no good if 
the figure-of-merit decreases. An iteration is considered better only if the 
figure-of-merit increases by a certain percentage (as specified by the user) . 

The user may wish to "reward" a certain performance measure only up to 
a certain desired value. In other words, up to a desired value, the figure-of- 
merit will include the actual value of the measure. When the measure exceeds 



this desired value, the f igure-of-merit will only include this desired value. 

If this happens for all performance measures, the figure-of-merit will not 
change at all from one iteration to the next, and the algorithm will have 
found an optimum that not only meets the requirements, but also satisfies the 
desired objectives. Examples 1 and 2 of Appendix B will illustrate cases where 
this happened. 

The case most likely to be encountered is when not all performance 
measures exceed their desired values, and the f igure-of-merit improves only 
slightly from one iteration to the next. The percent improvement is less than 
that required by the user and the design process terminates. For this case, the 
local optimum may yield a solution that satisfies all the design requirements , 
in which case the problem is solved. However, if the local optimum does not 
satisfy the requirements, the user must then either (1) try another initial 
condition for the autopilot variables, (2) add more complexity to the autopilot, 
and/or (3) relax some of the design requirements and/or alter some of the 
design objectives. 

The second way that the design may be terminated is directly by the 
user. He may specify termination after a certain number of major iterations or 
after a certain amount of computer time. 

2.11 Overall Flow Chart 

Figure 2-13 is a detailed overall flow chart of this nonlinear programming 
algorithm as it is applied to the problem of launch vehicle autopilot design. 

It shows the flow of information from the initial autopilot, all the way through 
the step-size optimization routines, through the termination routine, and back 
to the initial autopilot for the next major iteration. 

In the block showing where the partial derivatives are calculated, 


reference is made to closed-loop roots. This illustrates the possibility of 




Figure 2-13 














not only optimizing stability piargins, but also of optimizing locations of 
at least the so-called dominant closed-loop roots. In addition, a tolerance 
constraints routine is indicated. This routine is utilized to prevent the 
attitude error vector, the rate vectors, the accelerometer loop vector, etc, 
from getting larger than the total resultant vector at all frequencies. When 
the individual vectors become much larger than the resultant, vector cancellation 
may result, which can lead to problems when the airframe parameter tolerances 
are considered. 

In the block showing the constraint matrix set up, reference is made to 
the "Drift Minimutfi" condition, (Greensite £l4]and Hoelker||l7j ) which basically is 
a steady-state relationship between the autopilot gains. This steady-state 
relationship results in a mix between the forces due Co gravity, aerodynamics, 
and control deflections, that yields a zero net force perpendicular to the 
vehicle’s velocity vector. This illustrates that it is possible not only to put 
constraints on margins and root locations and on min-max values of the autopilot 
variables, but also to constrain the autopilot variables to satisfy other con- 
ditions like the Drift Minimum condition. 

Figure 2-13 also shows a block labeled "Simplex Algorithm." This is in 
reference to Dantzig's method [9]for solving linear programming problems.* 

The block labeled "Performance Objectives” refers to the fact that 
design objectives are used to form the cost function, while design requirements 
are used to form the constraint matrix. 

2.12 Sten-Size Ootimization Flow Charts 


Figure 2-14 illustrates the step-size optimization routines (namely, 
the Minor Loop and the Inner Loop) and their relationship to the other routines 
shown in Figure 2-13, 


*See Appendix E 



The Inner Loop is used to relax constraints in order to yield a feasible 
solution to the linear programming problem. The Minor Loop is used to 
maximize the autopilot variable step-size in order to get the "maximum mileage 
out of each set of partial derivatives, and in order to allow the algorithm 
to steadily converge to a local optimum, particularly to an interior optimum. 
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CHAPTER 3 


SUMMARY 


3„1 Suimnary and Conclusions 

As shown in the examples of Appendix B, the COEBRA program clearly 
demonstrates that this algorithm successfully solves the problem of auto- 
mating practical launch vehicle autopilot design and optimization. Perhaps 
the primary reason for the success of this algorithm is that its approach to 

design is much the same as the engineer's approach. 

Via this algorithm, the COEBRA program satisfies the five basic design 

requirements : 

(1-a) The COEBRA program deals directly with the "optimization” of 

performance measures in both the frequency domain and the time 
dom* in as well as to constrain the location of the rigid body 
dominant roots. COEBRA imposes an inequality constraint in 
each individual performance measure and each pair of dominant 
poles at all flight conditions being considered simultaneously. 
With the optimum performance measure cost function, COEBRA not 
only realizes the minimum performance measure requirement, but 
also seeks to "optimize" all performance measures. With a cost 
function distinct from the constraint matrix, the cost function 
can be formed from the performance measure Qb1ectives > while the 
constraint matrix can be formed from the performance measure 
requirements ; 

(1-b) With the "optimize load reduction" cost function (formed via a 
time domain transient response routine), COEBRA seeks to mini- 
mize structural bending moment loads ( /3 , S^j andS^ due to winds) 
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while meeting the minimum performance measure requirements. 

(1-c) COEBRA can constrain the autopilot parameters to the so-called 
Drift Minimum condition [l^ [l 7] , thereby minimizing trajectory 
dispersions. In fact, COEBRA can design a Drift Minimum auto- 
pilot that has the maximum amount of load relief capability 
and that meets the minimum stability margin requirements; 

(2) COEBRA designs with a user-selected autopilot configuration. 

From the outset, only practical controllers are considered since 
the user selects the number and types of feedback loops and the 
number of gair« and filters. COEBRA optimizes the values of the 
parameters within this feedback structure and constrains the 
minimum and maximum allowed values on each parameter; 

(3) COEBRA handles the problem of multiple time point design by 
forming the cost function and matrix of constraint equations 
from performance measures and wind responses at several time 
points or vehicle states for both deterministic and random forcing 
functions. In this manner, all vehicle states are optimized 
simultaneously. Autopilot parameters can be shared between the 


i: I 


vehicle states . 

A novel feature of the COEBRA program and this design algorithm 
is that multiple time points are handled by considering a 
separate airframe for each time point. It is obvious that these 
"separate airframes" can come from the same flight condition or 
trajectorv time point. For example, the "first airframe can be 
the nominal airframe at time t^^ while the "second airframe can 
be the airframe at time t^ with a tolerance on one or more of the 
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, COEBRA. can treat both the 


vehicle parameters. In this way 
nominal and the toleranced airframe together to yield a 
single autopilot that vjill handle both conditions; 

(4) COEBRA can handle a very high order system (30th and greater 
with up to eight bending and slosh modes per time point). 

With a user selected feedback configuration, the complexity 
of the autopilot does not necessarily increase with an increase 
in the order of the fixed parts of the system. Parameters like 
sensor and actuator dynamics are included in a very straight- 
forward manner and their inclusion only increases the required 

computations ; 

(5) COEBRA designs analog autopilots via the S-plane frequency 
response, and digital autopilots via the W-plane frequency 
response. 

Examples in Appendix B show that this algorithm can handle both 
interior and exterior optima. The examples also show that the 
initial conditions on the controller parameters need not yield 
feasible solutions, i.e., solutions that meet the constraint 
requirements. In fact, the examples of this appendix demonstrate 
that the initial condition on the autopilot parameters need not 
even yield a stable system. 

3^2 Projected Applications 

Even though the class of problems this algorithm can handle has not been 
established, it would appear that it can handle a large variety of engineering- 

type problems. 

For example, it would appear it can handle the problem of designing an 
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airplane flight control system with the so-called flying qualities design 
criteria[2] . These criteria include: (1) the longitudinal plane requirements 

on phugoid stability, flight path stability, and short period response; (2) 
the lateral-directional flying qualities criteria on the responses of the 
dutch-roll mode, the spiral mode and the roll mode; and (3) miscellaneous 
requirements on capability to perform crosswind landings, coordinated turns, 
etc. These criteria could simply be added to the flexible-body stability 
margin design requirements that are already included in COEiBRA, 

Another problem that this algorithm can handle is the design of a 
reaction control system. This type of control system uses discrete control. 

This algorithm can be used to optimize phase plane switching logic like the so- 
called "near-minimum-fuel" switching logic developed by Carney and Conover[5], 

Their phase plane logic was developed for a digital attitude control system 
that requires no rate gyros. 

Another problem that this algorithm can surely handle is the design 
of the autopilot for an interplanetary spacecraft like the Mariner [Kopf, 22], 

In order to handle the Mariner autopilot design problem, COEBRA uses a transient 
response routine that puts requirements such as rise time, overshoot, settling time 
and steady- state error on the vehicle's attitude response due to guidance 
commands , 

To conclude this discussion, a paper written by Robinson[28] is noted. 

In this paper, Robinson states that the COEBRA algorithm should prove fruitful 
in the optimal control of distributed parameter systems. 
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APPENDIX A 


COMPUTATIONAL AIDS FOR PAPER AND PENCIL DESIGN 


pbbcbding 




This appendix provides a few simple computational aids for a preliminary 
paper and pencil design, A table of rectangular to polar coordinate trans- 
formations (p A-2 to A-4 greatly facilitates the construction of the phase/ 
frequency plots, from which, with the corrected asymptotic gain- frequency 
plots (Bode plots) obtained with the assistance of the standard frequency 
response curves of pages A-5 through A-i6, the gain-phase (Nyquist) diagrams 
can be constructed. Representative gain phase curves for-«c<.co^ are 

shown on pages A-)7 through A-20. 

For completeness, pp. A-2l through A- 24 illustrate representative root 


loci for closed-loop systems. 
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APPENDIX B 


SAMPLE APPLICATIONS OF COEBRA 


This appendix presents some of the results obtained using the COEBRA 


program. The first two examples are classical testbook examples and the last 


five demonstrate the application of COEBRA to booster autopilot design. 


Throughout this appendix the following shorthand notation, typical of 


COEBRA formatting, will be used to represent system transfer functions in 


either the s-plane or the w-plane. The transfer function given by 


K (1 + T S) (1 + S + ) 

'• “l 


(1 + T S) (1 + T S) (1 + s + ^ 


will be abbreviated as follows: 


K (Tl, tl, wl) __ 
(T2, T3, r2,w2) 


Example 1. 


This ecample was selected to illustrate how COEBRA can be utilized to 


solve the digital system problem. Figure B-1 ic the block diagram of the 


system, together with performance measures to be achieved. The first step is 


to convert the z-fransforms to the w-plane. Thus, 


G, (w) = 0.13K (z + 1.31) (2 + 0.0_^ 

^ z (z-1) (z-0.368) 


1 + w „ ^ 
z = - , K«3 

1 - w 


= 0.77 g-w) g-0.13^w) (1-H).8925w: 
W(1 + w) (1 + 2.165 w) 


where K«3 is determined from the value theorem 


lim Gjjq G^(uT) = lim (z-1 ) Gj^ (z) 1,5 « • 

U — > CO z 1 




Digital Controller 


T * 0.5 sec, 

K >1.5 

y 

> 50 degrees 

M < 1.3 
P 


Figure B-1 Sampled Data System, Example 9-3, p. 291, Analysis and Synthesis of 
Sampled Data Control Systems, B, C, Kuu, Frentice-Hall, 1963. 










The digital controller configuration selected was of the form 

D(w) = 1 -h a 7 w ^ 

X-t-Tw 

Conventional paper and pencil design procedures, outlined in Figure B-2 yield 


a D(w) as 

D(w) = 1 + 50w 

1 + 100 w 


and an M^of 1»2, M of 50 degrees, of 12 db» converting .o the z-domain 


D(w) = 1 -t- 50 w 


1 + 100 w 


z = 


2-1 

2+1 


0,25 (z - 0,96) 
z-0.99 


0 • c 


and in the s domain, the compensation network corresponding to D(z) is 
evaluated from 

Z Gc(c-) D(z) 

s 


or 


Gc(s) = 


1 + 12.5s 
1 + 50 s 
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Figure B-3 illustrates the w-plane problem as submitted to COEBRA, 

Figure B-4 compares the results obtained from COEBRA with that obtained using 
a '".lassical pencil and paper approach. for the first COEBRA run* For this 
COEBRA run, both T1 and T2 were initialized with a value of 20, and hence 
the initial response shown in Figure B-4 is identical to the response of the 
plant onlyo Note that was 0,542 for the initial response. Figure £-4 shows 
the results of the first and sixth (final) major iterations of COEBRA. The sixth 
iteration was the final once since COEBRA was not "rewarded" for doing better than 
the classical design. In other words, recalling the discussion on termination, 
the figure-of-merit was not allowed to increase once the classical design results 
were matched. The results of the fifth and sixth iterations were identical, since 
it took COEBRA one iteration to decide that improvement was no longer possible or 
"permitted o" 

Table B-1 summarizes the initial and final compensators, as well as the 
final stability margins obtained from COEBRA run #1. Note that as with the 
classical approach, the final answer from Run 1 was a phase-lag compensator. 

Due to the circumstances of the problem, this minimum complexity (first order) 
compensator had to be a phase-lag model. In other words, phase- lead compensation 
would be ineffective. 

Additional COEBRA runs were made in an attempt to "map the hill," or in 
other words, to see what COEBRA would do with different initial compensators. 

As shown in Table B-1, Runs 2, 3, and 4 achieved essentially the same results as 
were obtained from Run 1 and the classical approach. With the initial compensator 
of FwUn 5, COEBRA climbed a local optimum that did not satisfy the design require- 
ments. Run 5 automatically terminated after 14 major iterations when the margin 
counter and the figure-of-merit essentially ceased to increase. Note that the 
final answer from Run 5 was not a phase-lag compensator. Run 6 was made with the 
denominator time constant (T2) of Run 5 changed to a value of 4., so that the 
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Figure B-4 Example #1, Gain- Phase Frequency Response 
Plot Comparing the Classical Design with 
COEBRA Run #1 
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inicial compensator was a phase-lag compensatoro Table B-1 shows that Run 6 
achieved essentially the same result as obtained using the classical approach^ 
Since the first order phase- lag compensator is the minimum-complexity 
compensator that can solve this problem, it is not difficult to understand why 
COEBRA could not converge to a final solution from every initial conditiuno 
This points out that the difficulty of any problem is dictated more by Hie 
degrees of freedom in the compensator than by the complexity or order of the 
plant, 

B-2 Example #2 

This example, easily carried out by classical methods as shown in 
Chapter 1, is included here to illustrate the use of COEBRA when the compensation 
is of higher order chan that of example 1. Figure B-5 shows the system in COEBRA 
format. The details of the classical approach are given in Chapter 1. The 
problem given to the COEBRA program was tc> adjust the four compensator time 
constants until the above three design requirements were satisfied. 

Table B-3 compares the classical results with those obtained from the 
first COEBRA run. With all the time constants in the compensator initialized 
to unity, COEBRA, in six major iterations, climbed to a local optimum that did 
not meet the design requirements. Conventional design procedures showed that 
the minimum-complexity compensator that is required to solve this problem, is a 
lag-lead compensator, and the final answer for the unsuccessful Run 1 is r|£t a 
lag-lead compensator, COEBRA terminated after six iterations when the margin 
counter and the figure-of-merit ceased to improve, 

COEBRA was reinitialized to the compensator shown for Run 2 in Table B-3. 
As can be seen from the table, Run 2 achieved all the design objectives. It did 
so with a lag- lead compensator. 

Figure B-6 compares the classical results with those obtained from Run 2, 
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Figure B-6 Example #2, Gain- Phase Frequency Response 

Plot Comparing Classical Results with COEBRA, 
RRun #2 
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Since Che numerator and denominator of the initial compensator for Run 2 are 
identical, Che inital response shown in Figure B-6 is identical to the response 
of the plant only. As can be seen, the system with a unity compensator (COEBRA's 
ini! i compensator) is unstable. Note Chat for Che initial lesponse is 22.4 
rad/scc. Figure shows the results of the fourth and sixth minor iterations 

in the first major iteration of Run 2, The results of the sixth minor iteration 
in the first major iteration satisfy all the design requirements, and these 
results were the best COEBRA was "allowed" to achieve. As with Example #1, the 
reason for this was that COEBRA was only "rewarded" up to the design requirements. 
In other words, recalling the discussion on Termination, the f igure-of-merit was 
not allowed to increase once the classical results were matched. Run 2 ran for 
two major iterations since it took COEBRA one iteration to decide that further 
improvement was not allowed. 

Two more COEBRA runs were made on this problem. Their initial compensators 
were chosen to "lie between" the initial compensators of Run land Run 2, From 
Table B-3, it is seen that the initial compensator for Run 3 was "close" to that 
for Run 1, As with Run 1, Run 3 climbed a local optimum that did not yield a 
feasible solution. Run 3 went four major iterations before the counter and 
f igure-of-mcric indicated that no further improvement was possible. The final 
answer from Run 3 was not a lag- leg compensator. 

Run 4 was made with an initial compensator that was "between" the 
initial compensators for Run 2 and 3, As can be seen from Table B-3, Run 4 
achieved satisfactory results with a lag-lead compensator. 

Since a second order lag- lead compensator is the minimum-complexity 
compensator that can solve this problem, it is not difficult to understand why 
COEBRA could not converge to a final solution from every initial condition. As 
with Example #1, this points out again that the difficulty of any problem is 
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dictated more by the order oi the compensator than by the order of the plant. 

In most cases ) the **optimurQ hill** boradens and smoothes out as the order or 
complexity of the compensator increases. 

Table B-4 summarizes the computer time requir-'' to nake COEBRA runs 1 
through 4. 

B-3 Example #3 

Example #3 illustrates the application of COEBRA to a single-time-point 
autopilot design problem where the initial autopilot was so poor that it resulted 
in a rigid-body instability. The objective of the COEBRA run was not only to 
stabilize the system, but also to optimize all stability margins. 

The airframe (or system to be controlled) included rigid-body dynamics 
and eight structural bending modes. The block diagram of the airframe/autopilot 
system is shov n in Fi.-jre B-7. Since this is a so-called analog autopilot, the 

design is performed in the S-plane. 

Figure B-7 shows the attitude loop with a gain and two filters, and 
two rate loops (for bending mode stability [^Harris, 15J), each with a gain and 
two filters. The problem given to COEBRA was to adjust these nine autopilot 
parameters until all stability margins were optimized. 

Figure B-8 shows the open-loop frequency response resulting from the 
initial autopilot. This figure shows that the initial autopilot did result in 
a rigid-body instability. The resonances of the eigh: rt.ructural bending modes 
are indicated by Figure B-8, which also shows that the rigid-body phase margin 

frequency (^^^) was 4.62 rad. per second. 

Arrows around the critical point in Figure B-8 illufc'trate the required 
rigid-body and first mode stability margins. It was also required that be 
greater than 2, rad/sec and that modes 2 through 8 be gain stabilized with their 
peaks resonating below "-10** decibels. 
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Figure B-7 System Block Diagram for Example #3 








Figure B-9 shows the frequency response after the first major 
iteration!, The system is now stable, with = 2,07 rad/sec. 

Figure B-10 shows the frequency response that resulted from the third 
and final iteration, GOEBRA self- terminated after all design requirements were 
met, and after the margin counter and f igure-of-merit ceased to significantly 
improve following the second major iteration. In other words, tne results of the 
second .md third iterations were identical since it took COEBRA one iteration to 
determine that design improvement was no longer possible. Further improvement 
would have been rewarded, but COEBRA was unable to achieve it. The margins 
that prohibited further improvement were and the phase margin on the "backside' 

of the first mode. 

Table B-5 summarizes the results of this example. It shows the valu.is of 
all nine parameters for both the initial and the final autopilots. Table B-5 
shows that a satisfactory design was achieved in 493 seconds of computer time. 

The following is a discussion of how COEBRA presently treats the require- 
ment on the dominant rotational rigid-body closed-loop roots. Up to the present, 
the time domain response due to guidance commands of a large aerodynamically 
unstable flexible launch vehicle, has not been too critical. The major concern 
has been with stability under tolerances and with structural bending moment loads. 
The main reasons for specifying dominant closed-lcop root locations have been to 
(l)keep the autopilot frequencies sufficiently separated from the guidance loop 
frequencies for stability purposes, and (2) merely provide "somewhat adequate- 
response to guidance commands. Historically, it has been found that if the 
rigid-body "hase margin is greater than a certain value, then the rigid-body 
rotational closed-loop roots will be sufficiently damped at a high enough frequency 
For example, on launch vehicles like the one represented in this example problem, 
if the phase margin is greater than 30 degrees, with a frequency greater than 2.0 



, Gain-Phase Frequency Response 
from Initial Autopilot (u) = 4. 


Example #3 
Resulting 



Amplitude Ratio (Decibels) 



Phase (Degrees) 


Figure B-9 Example #3, Gain-Phase Frequency Response Plot 


Resulting from First Major Iteration (u» = 2.07 rps) 
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Table B-5o Example #3 Sunuiary of Results 


Initial Autopilot 


Final Autopilot 


.0333 ,0333 


,018 ,018 


T3 

T4 

1CR2 

T5 

T6 

o2 

.8 

,2 

,2 

,2 

.108 

,22 

.108 

.108 

.108 


Objective: Maximum Margins 

Single Time Point, 27th Order System 

3 Major Iterations and 19 Minor Iterations 

Computer Time: 493 sec or 164 sec per Major Iteration 

Note: Initial Autopilot Yielded Rigid- Body Instability 






























radians per second, then it is almost certain that the rotational closed-loop 
roots will have a frequency greater than 1,5 radians per second and a damping 
ratio greater than 0,30, The M circles for unity feedback systems tend to 
indicate they this has been. Hence, rather than finding the actual roots, 

COEBRA treats the requirement on the rotational closed-loop roots by putting 
minimum allowed values on the rigid-body phase margin and its frequency 
This approach was taken in order to avoid the computer time required to find the 
actual closed- loop roots. 

It is recognized that the rigid-body response of a launch vehicle is 
comprised of a so-called first-order drift root as well as the second-order 
rotational roots [Creensite, 14, and Harris, 15], Hence, since the rigid-body 
response is third-order, the location of the rotational roots alone is not 
sufficient to ensure adequate response to guidance commands. However, since 
most launch vehicles are aerodynamica lly unstable, the instability of the vehicle 
generall serves to keep the attitude gain high enough, and the flexibility of 
the vehicle generally serves to keep the rate gain low enough, so that the rota- 
tional roots dominate the drift roots. In this way, the location of the rotational 
roots themselves can be used to indicate response to guidance commands. For this 
example problem, where the final autopilot yielded a phase margin of 38. degrees 
at a frequency of 2.03 radians per second, the closed-loop rotational roots had 
an effective damping ratio of 0*68 and an undamped natural frequency of 2*9 
radians per second. This satisfied the design requirements, and a transient 
response showed that these roots dominated the drift root which had a time 
constant of 11.5 seconds. 


*Nichols Chart of Appendix A 



B-4 Example #4 

Example #4 illustrates the application of COEBRA to the same airframe 
that was used in Example #3, but this time the initial autopilot was so poor 
that it resulted in a first-mode instability. As with Example #3, the objective 
of the CQEBRA run was not only to stabilize the system, but also to optimize 
all stability margins. 

The block diagram of the airframe/autopilot system is the same as that 
of Example #3, and is sho\m in Figure B-7 . The problem given to COEBRA was to 
adjust the nine S-plane autopilot parameters until all stability margins were 
optimized „ 

Figure B-11 shows the open- loop frequency resp'nse resulting from the 
initial autopilot. This figure show'' that the initial autopilot did result in 
a first-mode instability. The resonances of the eight structural bending modes 
are indicated in Figure B-11, which also shows that the rigid-body phase margin 
frequency (co ) was 1,17 rad, per second. 

As with Example #3, the arrows around the critical point in Figure B-11 
illustrate the required rigid-body and first mode stability margins. It was 
also required that w be greater than 2,0 rad/sec, and that modes 2 through 8 
be gain stabilized with their peak resonating below "-10" decibels. Figure B-11 
shows that with the initial autopilot, the fourth and fifth modes exceed this 
requirement. 

Figure B-12 shows the frequency response after the second major iteration. 
The system is now stable, with = 1,3 rad/sec. 

Figure B-13 shows the frequency response that resulted from the fifth 
and final major iteration, COEBRA self-terminated after all design requirements 
were met, and after the margin counter and f igure-of-merit ceased to significantly 
improve following the fourth major iteration. In other words, the results of the 
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Figure B-lla Example #4, Gain-Phase Frequency Response Plot 
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Figure B-13. Example #4, Gain-Phase Frequency Response Plot 
Resulting From Fifth (Final) Major Iteration 
(cj - 2,0 rps) 










Table B-6o Example #4 Summary of Results 


KRl I T3 


T4 KR2 


Initial Autopilot 


Final Autopilot 


.0333 .0333 


.093 .061 


.0333 .0333 1 ,5 


,0333 .0333 


.053 .32 .039 .021 


Objective: Maximize Margins 

Single Time Point, Same 27Ch Order System as Example #3 
5 Major Iterations and 29 Minor Iterations 
Computer iime: 975 sec or 193 sec per Major Iteration 

Note: Initial Autopilot Yielded First Mode Instability 
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fourth and fifth iterations were identical since it took COEBRA one iteration 
to determine that design improvement was no longer possible. Further improve- 
ment would have been rewarded, but COEBRA. was unable to achieve it. The 
margins that prohibited further improvement were , the phase margin on 
the "backside" of the first mode, and the modal peaks of the third and fifth 
modes. 

Table B-6 summarizes the results of this example. It shows the values 
of all nine parameters for both the initial and the final autopilots. Comparing 
these values with those obtained from Example #3, it is seen that COEBRA "climbed" 
to a different local optimum than it did for Example #3, However, some features 
of the two results are similar. For both, the attitude gain (KD) is about 1,5, 
the total rate gain (KRl + KR2) is about unity, and KRl is greater than KR2 
in order to "center" the first mode around zero degrees. 

Table B-6 shows that a satisfactory design was achieved in 975 seconds 
of computer time, 

B-5 Example #5 

Example #5 illustrates the application of COEBRA to a three-time-point 
autopilot design problem, where the objective was to optimize structural bending 
moment load relief capability, COEBRA was initialized with an autopilot that 
had previously been designed by engineers. The reason for tiiis COEBRA run was 
to determine if design improvement could be achieved. Design improvement was 
defined as an autopilot that had more load relief capability, but still met 
satisfactory stability margins. 

The example deals with the max-q portion of flight where aerodynamic 
loads are critical. The three vehicle states that are designed together are: 

(1) the time at which the lead relief loop (the accelerometer feedback loop) 
is switched in; (2) the max-q time point; and (3) the time at which the load 
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relief loop is switched out. The airframe consists of rigid-body dynamics 
plus three bending modes at each time point. 

Figure B-U shows the block diagram of the airframe/autopilot system. 

In addition to the attitude loop and two rate loops. Figure B-14 shows the so- 
called load relief loop. This is a feedback loop on a lateral body-mounted 
accelerometer signal, and Harris [iS] discusses how this loop is used to reduce 
the angle of attack and control deflections (hence bending moment loads) in the 
presence of the wind, 

COEBRA. is allowed to vary the gains and filters sho\m ir Figure B-14. 

This is a digital autopilot design problem, and hence, these gains and filters 
are defined in the W-plane. Of course, when the design is complete, these gains 
and filters will be transformed to the Z-plane where they will be mechanized as 

coefficients in different equations. 

Figure B-14 shows that, at each time point, 15 autopilot parameters can 
be varied. Since this is a digital autopilot, the fcur gains (KD, KL, KRl and 
KR2) can be different at each of the three time points. The 11 filter network 
values, though they may be varied, must have the same values at all three time 
points. 

Figure B-15 shows the open-loop frequency response plot that results at 
the max-q time point from the engineer's final autopilot (initial autopilot for 
the COEBRA- run). This figure, as well as the frequency responses at the other 
two time points (not shown), show that all margin requirements are satisifed. 
When a six-degree-of-£reedom (6 DOF) trajectory was run using this final 
autopilot, the load relief indicator (which is a product of the dynamic pressure 
times the angle of attack and is indicated as q a) was 4908 pounds per square 
foot. 

Figure B-16 shows the frequency response plot that resulted at the max-q 
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Example #5, Gain-Phase Frequency Response PloC_ 
Resulting from Third (Final) Iteration at Max-q 
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Table B-7c Example #5 Summary of Results 


Objective: Maximize Load Relief Capability 


Autopilot 

Stability Margins 

qa 

COEBRA Initial 

Satisfactory 

4908 

(Engineer's Final) 



COEBRA Final 

Satisfactory 

4765 


System Order: e 28th Order at Each of the Three 

Time Points 

o 15 Autopilot Variables at Each 
Time Point (4 gains can have 
different values at each, 11 
filters must have same values 
at each). 

3 Major Iterations 

Computer Time: 21„3 Minutes or 7-1 Minuters per 

Iteration 
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time point from the third and final iteration of the COEBRA run. Stability 


margin requirements ar met at this time point, as well as at the other two 
time points. For this final COEBRA autopilot, a 6 DOF trajectory simulation 
showed that q had been reduced to 4765 pounds per square foot. 

Figure B-16 shows that the stability margins from COEBRA’s final auto- 
pilot, though satisfactory, are less than those from the engineer's final autopilot 
(Figure B-15). This demonstrates the tradeoff that does exist between stability 
and load relief. 

The conclusion of this example is that, starting from the engineer's 
final autopilot, COEBRA was able to achieve an improved design by adjusting the 
values of gains and filters within an engineer's established configuration. It 
is noted that t le COEBRA improvement in load relief did not result because 
the engineer was incapable, but rather because he was not required to obtain more 

load reduction. 

Table B-7 summarizes the results obtained from this example. It shows 
that a satisfactory result was obtained after 21.3 minutes of computer time, or 
7,1 minutes per iteration. 

One final note is mentioned at this time. Another COEBRA run was made 
on this problem, beginning with the same initial autopilot, but with the objective 
changed to maximizing stability margins instead of load relief capability. The 
result of this second run was a design with improved stability margins, but with 
reduced load relief capability. Computer t.me for this run was 2.3 minutes per 

iteration, 

B'-6 Example #6 

Example #5 demonstrated the effectiveness of the structural load relief 
optimization phase of COEBRA when the initial autopilot met all of the margin 
requirements. Example #6 was run to see how the load relief phase performs when 
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the initial autopilot does not meet the margin requirements. 

The initial autopilot for this COEBRA. run was obtained as follows „ 

In several booster autopilots, there is a feedback loop that is used soley for 
high frequency stabilization. This loop is “washed out" at frequencies below 
the rigid-body phase margin, and serves to compensate for the load relief loop 
gain at high frequencies so that the load relief loop gain can be increase. So 
for Example #6, a previously designed autopilot case was chosen, and this "high 
frequency" loop was zeroed out. This resulted in unacceptable stability margins. 
The design objective for the COEBRA. run was to not only return to the condition 
where all margins are met, but also to achieve at least the same amount of load 
relief that was achieved with the engineer's original autopilot that used this 
high frequency feedback loop. 

For example #6, three flight conditions were designed together: (1) load 

relief switch-in; (2) max-q; and (3) load relief loop switch-out. The airframe 
included four bending modes at two of the time points, and three modes at the 
third. 

Figure B-17 is a block diagram of the airframe/autopilot system. This 
is an analog autopilot design problem, and therefore the design is performed in 
the S-plane. There are 14 autopilot variables, but since this is an analog 
autopilot, each of these variables must have the same value at all three time 

points , 

Figure B-18 is the open-loop frequency response plot that resulted from 
COEBRA.' s initial autopilot at the max-q flight condition. It shows that not all 
margin requirements are met. This same situation exists at the other two flight 
conditions (not shown). 

The load relief indicator for the engineer's original autopilot that 
used the so-called high frequency feedback loop was 4490 pounds per square foot. 
This result was obtained from a 6 DOF trajectory simulation. COEBRA.'.s first step 
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was to "get feasible," but in so doing, it had to give up load relief 
capability. COEBRA. met all margins after three iterations, but qa (from a 
6 DOF simulation) increased to 4580 pounds per square foot. However, from 
the 4th to the 8th Iteration, all margin requirements remained satisfied, and qa 
began decreasing, until on the 8th and final iteration, it had decreased to 3975 
pounds per square foot. Again, this was obtained from r 6 DOF simulation, and 
this qa was 127. less than that of the original autopilot with the high frequency 
feedback loop. Figure B-19 does show that all margins are met at the max q 
time point with the results of the 8th and final COEBRA iteration. The same 
situation existed at the other two flight times. 

The folloxi;ing discussion refers to Section2.9 of Chapter2 on Convergence 
to an Exterior Optimum. This example has dramatically demonstrated how the 
COEBRA algorithm converges to a constrained optimum from an unfeasible initial 
point. The first three iterations were required in order to reach a feasible 
solution. In "getting feasible," load relief capability was reduced. This did 
not necessarily have to happen, since the algorithm does try to optimize while 
"getting feasible." Once the feasible region was reached, the algorithm moved 
along or parallel to the constraint boundaries until the constrained optimum 
xvas reached. The fact that this actually occurred is known because load relief 
capability steadily increased from the 4th to the 8th iteration, while several 
stability margins remained "tight against" their requirements. Two of these 
"tight margins" can be seen in Figure B-19, These two margins are called the 
rigid-body phase margin, and the phase margin on the "backside" of the first 
structural bending mode. These margins are indicated by arrows in Figure B-19, 
These margins xi/ere tight after three iterations, and remained tight from the 4th 
to the 8th and last iteration. 

Table B-8 summarizes the results of Example #6. This table shows that 
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Figure B-19, Example #6, Gain-Phase Frequency Response Plot _ 
Resulting From Eighth (Final) Iteration at Max-q 




Objective: Maximize Load Relief Capability 


Autopilot 

Stability Margins 

qa 

Engineer's Final 

Satisfactory 

4490 

COEBRA Initial 

Unsatisfactory 

4490 

Third Iteration 

Satisfactory 

4580 

Eighth (Final) 
Iteration 

Satisfactory 

3975 


• System Order: o 25th Order at 2 Time Points 

o 23rd Order at 1 Time Point 

o 14 Autopilot Variables Which 
Must Have Same Value at Each 
Time Point 

• 8 Major Iterations 

• Computer Time: 59.2 Minutes or 7.4 Minutes per 


Iteration 








the computer time for this example was 59«2 minutes, or 7=4 minutes per 
iteration* 

B-7 Example #7 

Example #7 illustrates using COEBRA. to design a load relief autopilot 
in two phases: (1) the initial phase being to first meet all margin require- 

ments; (2) the second phase being to optimize load relief capability* For this 
example, this approach was considered essential because the "first guess" or 
initial autopilot was very poor. 

This example is taken from a recent effort to design an autopilot for a 
space shuttle booster configuration. COEBRA was used to design the autopilot 
for all three channels (pitch, yaw, and roll) at all the critical flight 
conditions during the first two minutes of ascent. At all the flight times, 
the airframe included from seven to eight structural bending modes. 

The flight condition for this example is the yaw channel during the max' 
portion of flight. Three time points were designed together: (1) load relief 

switch-in; (2) max-q; and (3) load relief switch-out. The airframe included 
seven modes at each of the time points. While all the results obtained from 
this design effort are worth noting, this example was selected since it illus- 
trates the two phased approach to load relief autopilot design. 

Figure B-20 is the airfrarae/autopilot block diagram for this example. 

It shows the attitude loop, a rate loop, and the load relief loop. In addition, 
it shows an attitude acceleration loop. This is the so-called high frequency 
loop that was referred to in Example #6, Figure B-18 shows that there are 20 
autopilot variables at each time point. The four gains can have differ -nt 
values at each time point, but the 16 filter parameters must have the same value 

for all the time points. 

Figures B-21, 22 and B-23 show the frequency response plots for the 
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Figure B-20 System Block Diagram for Example #7 















-40 


-60 


220 260 300 340 20 60 100 140 180 220 260 

Phase (Degrees) 


Figure B-21, Example #7, Gain- Phase Frequency Response Plot 
Resulting From Initial Autopilot at Load Relief 
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Figure B-23. Example #7, Gain-Phase Frequency Response Plot 
Resulting From Initial Autopilot at Load Relief 
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initial autopiloto The system is stable, but the initial autopilot is very 
poor„ The basic margin requirement is that modes 3 through 7 be gain 
stabilized with a peak amplitude below "-10" decibels. Only the first and 
second modes can be phase stabilized, but if they are gain stabilized, their 
so-called "closest approach" distance to the "-1" point must be equivalent to 
10 decibels. 

The first COEBRA run was made to optimize stability margins. After one 
iteration, all margins were met. The next COEBRA runs were made to optimize 
load relief. After six more iterations, an autopilot resulted that yielded the 
plots of Figure B-24, 25 and B-26. All margin requirements are met. A 6 DOF 
trajectory simulation was not made, but estimates based rn linear transient 
response results indicate that bending moment loads were reduced 25% from the 
initial to the final autopilot. Computer time required to do this job was 98 

minutes, or 14 minutes per iteratir'n. 

This example points to another way in which the COEBRA algorithm can be 
used. By observing the progress it is able to make from iteration to iteration, 
it can be used to design a minimum-complexity autopilot. For this example, the 
fact that COEBRA was able to satisfy all margin requirements in only one iteration, 
tends to indicate that some of the degrees of freedom in the autopilot could 
probably be eliminated. 

B-8 Conclusions 

The results presented in this appendix clearly demonstrate COEBRA' s 
ability to successfully design autopilots for large flexible launch vehicles. 
Experience with the program shows that while it generally does not save computer 
time, it does save manpower and the time required to design an autopilot. 
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APPENDIX C 


RANDOM FORCING FUNCTIONS 


In the foregoing analysis, though not explicitly stated, it has been 
assumed that the vehicle forcing functions were deterministic in nature. For 
example, the wind forcing function has been, in general, implemented by means 
of a pre-programmed wind profile based on a deterministic combination of steps 
and ramps. A more general characterization of the wind profile would be as a 
stochastic or random process. The vehicle behavior under the influence of such 
a forcing function can significantly affect the performance measures. Furthermore, 
many deterministic signals have random disturbances or noise superimposed on them. 

COEBRA has incorporated in its repertoire the capability to examine the 
effects of a random or noise forcing function. The prime occurrence in our case 
is that of a random wind profile. Consider, for example, the structural bending 
moment load, l(t), existing at a specific sensor location on the vehicle. This 
load can be expressed as 

1 (t) = K^a(t) + I . 

where is the vehicle angle of attack and 5^ is the control actuator function. 

The approach adopted by COEBRA is to constrain or optimize the mean square value 
of l(t) resulting from the random forcing function. 

Once again we call upon the o^en literature for a detailed exposition 
of the random process problem ^ The fundamental mathematical 

tools and the two-sided Laplace transform and the Fourier transform. Following 
the pattern of Newton, Gould & Kaiser , we make no formal distinction between 
them other than that s = jwin the foliating. Let the vehicle or system function 


be 


H(s) = 


C-2 


where f(t) is the random forcing function and H(s) thus contains the control 


C-1 


equations and vehicle dynamics or simply those system parameters to be 
considered or optimized. 

Associated with the system represented by H(s) and, of course, the 
random forcing function f(t) and the system output, l(t). The convolution 

integral yields 

/ CO 

h (t) f(t-T) dr 

Introducing the”concept of correlation functions, ensembles and the 

hypothesis as fundamental properties of systems subject to random forcing functions, 
we treat first the concept of an ensemble. To do so, we follow the influence of 
Newton, Gould and Kaisert^^to first make a subtle distinction between a random 
process and a stochastic process as used herein. A random process is one describable 
only in terms of its statistical properties. For our purpose, a stochastic 
process or signal is one that exhibits a degree of randomness but is not 
exclusively random; i.e., perhaps a sinusoid accompanied by noise. An "ensemble" 
is considered to be a set of stochastic signals, each generated by an Identical 
process. Viewed from the time domain point of view, the appearance of each 
element of an ensemble over the same time span would be different but the 
properties of the ensemble can be identified. For example, the average value 
of the product of a signal at a specified time, b^, multiplied by the value of 
the same signal at a different time, b^ +7:. is known as the ensemble average 

or auto correlation function, defined by 

(t ,r) t c-4 

If we denote the ensemble average of the stochastic signal , 1/'(t) as 




• • • ♦ 


whereA/r; is one element of the ensemble , and the time ..verage of the function 

V oO 


. Z-/ — ' 
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we observe that if 


the signalof^f) can be considered ergodic, a usual assumption when considering 
engineering systems. Using this ergodic hypothesis we can rewrite equation C-4 


as 


oo 




AMT 
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Returning now to our autopilot problem, defined by H(s), we note that 


oO ^ 

r 




» « 6 0 


C-8 


“ Otj ^ (Ui 

where is the autocorrelation function of the stochastic forcing function f(t). 
We now consider a frequency function related to the correlation function, 


— 0^ 


C-9 


Recall now, equation C-8, multiplying both sides by^ ^*^and integrating on jf* from 


to ^ we have 

cO 


f° 


C-10 


Now changing the order of integration on the right of C-IO so that we 
integrate with respect to T first, and adjusting the argument of the 


exponential to agree with that of the correlation function, we obtain 

ob eb 

- 






~ ci> 
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But this is just 
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C- 12 


Now, sincePJ^} given by C-9 is in truth the Fourier transform of 


we can write 


Z7T^ 


ja^ 

st 

$ ks)C Js 


-J cli 


• • • • 


C-13 


The integration is along the imaginary axis of the s-plane enabling us 'to 
write the inverse transform in terms of the real frequency or 



— c* 


Letting T= 0, we observe 

aO 


C-15 


But this is, from C-7, simply the mean-square value of l(t), the structural 
bending moment load# 

Thus, using C-12, C-15 and the knowledge of the frequency function 
we are able to make use of a constraint on the mean-square value, as the case 

arises, to assist in the design of the autopilot# 

A word or two about the frequency function ('■5; . Consider, for a 
moment, that f(t) is a stochastic voltage impressed on a 1 ohm resistor. The 
mean-square value of f(t) then represents the average power dissipated in the 

C-4 


resistor and 4* (j ^ ) represents the power-density spectrum of the signal f(t) 
Common practice refers to any frequency function obtained through the Fourier 
transform of an auto correlation function as a power-density spectrum. 

Using the COEBRA program, we can determine, as before, the gain 
constants and filter parameters such that l^(t), or a similar random response 
function is constrained and/or minimized. In other words, with <p specified by 
the user COEBRA will compute I^(t) and H(s) from the total dynamics. 



APPENDIX D 


INPUT FORMATS 

There are several fundaraental restrictions on the COEBRA design 
algorithm. If analysis is to be carried oat in the frequency domain with 
block diagram format, the autopilot configuration must meet that of 
Figure 1-7 (the QD030 subroutine). Within limits simple block diagram 
transformulations should lead to an equivalent structure of this form as 
illustrated in Figures D-1 and D-2. In the event that the initial confi- 
guration is such that the block diagram transformations appear to be 
unwieldly or tend to obscure the effects of certain autopilot parameters 
an alternate procedure has been provided. This procedure is closely 
related to the state variable description of a control system. However, 
one should not construe, in what follows, that state variable techniques 
are being applied. Referring to Figure D-1 the reader will note the 
identification of input and output signals, to the individual blocks 

or transfer functions. This identification is the first step in the 
formulation at the so-called matrix input description for COEBRA. Before 
going any further, we observe the restrictions on the elements of the 
matrix to be developed. Each element is to be limited to a second order 
polynomial ratio in S. The matrix for a single flight condition can have 

only one active forcing function. 

If analysis is to be carried out in the time domain, the maurix 
format must be used. This leads to a somewhat restrictive definition of 
a flight condition, as illustrated in Figure D-3. A new COEBRA input 
must be supplied for each flight condition. 
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FIGURE D-3 FLIGHT CONDITION CONCEPT 
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Examining Figure D-1 we can write the following set of equations 
involving only second order functions in s by the introduction of slack 
variables and formatted in the matrix notation us. 

Ax-B0 

where the vector B is comprised of a column of zeros exclusive of a single 

forcing function element. The solution of this equation yields the time 

response of the autopilot/vehicle system. 












In matrix notation this becomes 
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‘26,26 


= -1 
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‘27,23 “ 1 + _Sr 
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‘28,27 
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a ts — 1 

^28,28 


‘29,6 


^—T + 2 "T~ S + 1 


T 

-257 
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257 


257 


‘29,29 


= -1 


= 1 


It should be noted that this example involves a matrix solution for a matrix 
of order nxn where n is dependent on elements a^^ gj ^26 6’ ^29 6 involve 

the transformation T , the vehicle dynamics from input x^ to the individual 
station sensors, in this case at stations 88, 230 and 257 for example. 

A further caution must be observed. The slack variables ^ 2 ^, 

X are introduced since the COEBRA program was designed to handle only 
linear and quadratic terms in s. 


In the frequency domain the matrix for mat is somewhat more compli- 
cated but also as versatile. The format is fundamentally the same with the 
vehicle dynamics shown in a more explicit fashion. The usual procedure 
is to open the autopilot loop, typically at as shown in Figure D-4. 



^18 

Figure D-4. Opening the Autopilot Loop 


The individual vehicle parameters now become, one-at-a time, the system 
forcing functions. The equation 

AX = 


now appears as in equation D-4. 
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Where represents the coefficients of the vehicle dynamical equations, 

T the vehicle control equation coefficients for those control effectors not 
used as a forcing function and T^2 the coefficients of the autopilot equa- 
tions. the vehicle response variable resulting from the active forcing 

function, the control effectors not active as a forcing function, the 
autopilot variables as in equation D-2; b^^ is the three element vector as 
coefficients of the active control effector under study. vehicle 

submatrix will include 2 rigid body equations involving the lateral and 
rotational motions of the vehicle and n bending modes, thus is an (n+2) x 
(n+2) matrix; T. « is m x m where m is the number of control effector loops 

^ i. 

that remain closed (in active control effectors); least for the 

example of Figure D-1, 29 x 29 and is dependent not only on the actual 
autopilot configuration but also on the number of slack variables required. 

In formatting for input to COEBRA the ^ vector is used to augment 
the A matrix. If A is n x n, the column vector B is added to A as the 
n + 1st column as an input to COEBRA. The precise details are left to 
Volume II, the COEBRA Users Manual. 


appendix e 


THE SIMPLEX ALGORITHM 


The simplex algorithm is basically a simple procedure for solving a 
linear programming problem. The basis for this appendix is derived from 
references 49, and 50. 

The general linear programming problem can be stated as follows; 
determine the set of n nonnegative variables, x^, that satisfy m linear 
inequalities or equalities (the contraints). 


A X 

= m X n — n 


B , m n 
— m 


!- = -1 

where one and only one of the signs in the brace, = or holds for each 
individual constraint but can vary from one constraint to another within the 
set of m constraints defined by E-1. In addition, these n nonnegative varia 
bles are to extremize (maximize or minimize) a linear form (objective or 
cost function) 


z 



Any set x which satisfies E-1 is called a solution. Any nonnegative solu- 
tion set X ^ o is called feasible (admissable) and any such set that 
extremizes E-2 is called an optimum feasible (admissable) soulution. 

Generally, it is easier to work with equations than inequalities 
so we convert the inequalities through the introduction of slack variables. 
For example, given b^<o 


.il^i + a.2K2 +. . +^„^„j>=< 


multiplying through the inequality by -1 yields 


n 

-Z 


a . X . > - b , 


b. > o 






where multiplication b> -1 has reversed the inequality. For example, 
but (-7)(-l) ^ (“2)(“1)* is thus simple to convert* all constraints with 

each o. 

Converting the inequality constraints into equalities we observe the 
convention in ordering the constraints as follows: 

1. All constraints with ^ signs 

2. All constraints with signs 

3. All equality constraints 

Examining the type 1 constraints first we observe for the i-th constraint 


n ^ 

•? , a. . X. Ss. b. 

J==l iJ J 1 


Introducing a slack variable x^_j_^ _> o we have 


X 


n-f-k 


b. 

1 


n 

- L 
j=l 


a . . 

3-J 


X . 

J 
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Re-arranging E-6 we can write 


n 


ij J 


For a type 2 constraint 

X b 


n 

P=1 P P 


h 


We write 


n 

= a 

‘n+h P =1 h. 


X, 


b. 

r 



To obtain 


Thus, the original constraints are rewritten, with the introduction of slack 
variables, as 


a . . X . + X _ , 
ij J n+k 


E-Ua 


2^ a, X - X ,, 

, h,. n+h 

P=1 P P 


E-llb 


X a , X, = b, 

k=l ^ ^ 


E-llc 


It is thus possible to recast the linear progranming problem in the form 


A X = B 
=mxm — n —m 


E-12a 


x^ Oy (i = 1, 


E-12b 


max z = c X 


E-12C 


We assume in what follows that we have a basic feasible solution at hand and 
proceed to develop the algorithm to generate the "optimum" basic solution 
and then treat the problem of determining the first basic solution. 


= col =^20 - • • ■ ’'mo ° ° ^ 


be a basic feasible solution to the linear programming problem and the 
set of linearly independent basis vectors be P^, whence 

m 

£ X P.P. =P 

.1 10 “1 -1 -o 

J=1 


E-13 


E-14 


X, c. 
10 1 


E-3 


Where all > o, the are the coefficients from E-12c and the 
corresponding value of the objective or cost function. Since the are 
linearly independent any of the vectors ^2 • * • * > 
in terms of ‘ ‘ 


and define 

m 

z. = X X c , j = 1, . • • •» O 

^ i = 1 

where the are the coefficients corresponding to the we have - 

c.^o, evidencing a feasible solution. 

Denoting the set as 

i = •••••• 

we observe that Pq 

^ 5o = £0 


Inverting 



each of which is a column vector. 


We group the vectors in the form 


r P : P. : : P : P 4.1 • • • • 1 

L—o “1 *“2 ~m — nrrl n -i 


C In* ^ • -^+1' 


o '= —m+l 
Pre-multiplying by ^ ^ we- obtain 

-1 


: P 1 
— n -* 


S C£o= ^ 


: P 1 

“n -* 


= r X : I : Xj-n • • • ^ ~1* 

L — o mxm -Tn+1 ”tU 

Knowing the we form and examine for all j, to find z. - ^o. 

If, for all j, ]> o, we multiply E-16 by some constant 0 and 

subtract from E-14; similarly multiply E-17 and subtracting from E-15 we 
have 


m 

X (x. 

10 


i=l 


Ox. .) P. 4- dP. 
jj: “3. -J 


p 

— o 


and 


m 

X 

i*l 


(x. - ^x.,) c. + ^c. - z 

io 3-j 3- J ^ 


- 0 (z . - c .) 

J J 


where ^c. has been added to both sides of E-27 . If, in E-26, all the 
J 

^ 0 we have obtained a new feasible solution whose objec- 

^ io ij 

tive function value is, from E-27, z = ^ (z. - c.)^ j 


mVi ay- a 



If, in initially forming z. - c^, we find for some j, 

>je already have the minimum feasible solution. 

If the set of n vectors la’* * ‘ ^ contains m unit vectors 

that can be regrouped to form an m x m unit vector we write E-18 as 

^ 


Under this condition equation E-20 becomes 


X = P 
— o — o 


X. = P. 
-J “J 


and we have an initial solution 


To start the simplex procedure, in this case, we have x^^^ = b^, and 
To obtain the for j = 0, 1, . . . n, we form the inner 
product of the j-th vector with the given c, or 


z - E c. X. 

0 .1 1 10 
1=1 


= Z -=i ’'IJ- j = 1. 2, 
J 


It is convenient to arrange the elements of the problem in a tableau as 
shown in Figure E-1. 

Referring to Figure E-1 we enter the elements and z^ - c^ into 
the m + 1st row of their respective columns. If all the z^ - J “ 

1, 2 . . . n, the solution x^ = col ^20’ * ’ ' ^mo^l " ^^1’ ^2 

b I is a minimum feasible solution and the corresponding cost 
* ’ * ’ m— ' 

function is z . If not, we compute a new feasible solution whose basis 
o 

contains m - 1 vectors of the original basis P^, effi- 

cient procedure is to select a vector with its - c^ ^ o which yields 
the largest decrease in the value of z. This vector, should be the 

























one corresponding to the max 0^ (z^ - c^) where 0^ is defined by equation 

E-28. An alternate choice, often used, is to select that P. corresponding 

to max (z. - C.). This is the approach we take here. Thus, we let 
J J 

m^ (Zj - Cj) = - Cj^ 

whence the vector is to be introduced into the basis set. Computing 


X. 

Q = min x., 

° i 

where x., belongs to the x in the basic solution. For exanple 

Q = — 

° ^Ik 

and the vector P^ of Figure E-1 is the one to be eliminated from the basis 
set. The new feasible solution has as its basis P^, ' * * -1-1’ -1+1’ 

• * *=> -k* 

Since our initial basis is £ ^2 ' ’ * ^ T ^ 

we can express all of the vectors P^ in terms of this set, or 

= ""lo^l -1 ""nio 4 

4 = ^ik ^1 ’'.fk 4 ’'mk4 

Pj = X. . Pj^ + . . . + xi. P^ + . . . + X_^ 

Solving E-36b for P^ we have 
Substitution of this value in E-36a yields 

4 = =^io 4 + • • • + =^^0 ^4 - ^‘ik 4 ■ • • • 

D^ 1 a. ^ . 6 -ff+1 ^mo 
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+ (X - -^ X , ) P 
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which is 


equivalent to E-26 with j - k and B = - B 


The new feasible solution 
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By substituting E-37 into E-36c we obtain each not in the new 
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. m 
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basis 


set as 



where 

X. . = X, . - X., j 

i-J 3.J Xjjj^ ik 

X* ^ ^ 

n«k 

and Xfl, is denoted the pivot element. 

Now 

z' y X. . c. - 

I 

which, by substitution of E-40 for x in 
m 

z’ * X C. X. , 

a ^ ^ 1 xo 

1 i 

gives 

z' = z - ^ (z - c ) 

0 0 x^j^ k k' 

I • 

To obtain the new solution, x^, the new vectors _x^ and the corres- 
ponding z. - c. we transform every element in the tableau of 

^ ^ t i 

Figure E-1 by equation E-42 for and xjf. - . The new cost 

function is 

f ( 

^o " \+ 1,0 

and 

I 1 

Z. -C. = X - . 

j J in + 1, J 
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Once an initial computational tableau has been constructed, the 
simplex procedure calls for the successive application (l>ei, an iteration) 
of: 

1. The testing of the z. - c. elements to determine v/hether a 
minimum solution has been founds i*e., whether Zj - 

^ 0 for all j , 

2. The selection of the vector to be introduced into the basis 

if some z - c ^0, i.e., selection of the vector with 
J J 

maximum z . - c . . 

J J 

3. The selection of the vector to be eliminated from the basis 

to ensure feasibility of the new solution, i.e., the vector 
with min ( x^^^ for those x^^ ^ 0, where k corres- 

ponds to the vector selected in Step 2. If all x^^^ ^ 0, 
then the solution is unbounded. 

4. The transformation of the tableau by the complete elimina- 
tion procedure to obtain the new solution and associated 
elements . 

Each such iteration produces a new feasible solution and eventually 
yields the minimal or optimum solution or reveals an unbounded one. 

A simple exanq^le will serve to clarify the foregoing. We wish to 
minimize the cost function 


z = c 


^x=Q0 1- 302003 


subject to the constraints 
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b) STEP lA 
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In forming the set of the second step tableau of Figure E-4 we 


observe 
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IJ LJ X, 
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The are formed as before and the max (z . - c .) = z„ - c = 1/2 >0 


is observed. Thus ^ ^ ^ 
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= min 
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since we require x ^2 ^ 





is now the vector to be introduced into the basis set and the 


e to be eliminated to obtain the tableau of Figure E-5. 
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FIGURE E-5 STEP 3 IN THE COMPUTATIONAL SCHEME 


Since the max (z . - c.) = 0 this is an optimum solution. If the given 
linear programming problem does not contain a unit matrix, as was 
assumed in the foregoing the method of artificial basis detailed in 
references [.50 and 50 is referred to and will not be detailed here. 

It is the simple^ algorithm tliat COEBRA utilizes, as an effi- 
cient method of solving the autopilot problem* 

























